Counting Complex Disordered States by Efficient Pattern Matching: 
Chromatic Polynomials and Potts Partition Functions 
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Counting problems, determining the number of possible states of a large system under certain 
^ ' constraints, play an important role in many areas of science. They naturally arise for complex 

^— ■>) I disordered systems in physics and chemistry, in mathematical graph theory, and in computer science. 

. Counting problems, however, are among the hardest problems to access computationally. Here we 

suggest a novel method to access a benchmark counting problem, finding chromatic polynomials 
of graphs. We develop a vertex-oriented symbolic pattern matching algorithm that exploits the 
equivalence between the chromatic polynomial and the zero-temperature partition function of the 
Potts antiferromagnet on the same graph. Implementing this bottom-up algorithm using appropriate 
' computer algebra, the new method outperforms standard top-down methods by several orders of 

^SJ , magnitude, already for moderately sized graphs. As a first application we compute chromatic 

polynomials of samples of the simple cubic lattice, for the first time computationally accessing 
three-dimensional lattices of physical relevance. The method offers straightforward generalizations 
\ to several other counting problems. 

O . 

Given a set of different colors, in how many ways can one color the vertices of a graph such that no two adjacent 
vertices have the same color? The answer to this question is provided by the chromatic polynomial of a graph [1, 2], 
• which gives the number of possible colorings as a function of the number q of colors available. It is a polynomial in q of 
^ ' degree N , the number of vertices of the graph. The chromatic polynomial is closely related to other graph invariants 
: e.g. to the reliability and flow polynomials of a network or graph (functions that characterize its communication 
capabilities) and to the Tutte polynomial. These are of widespread interest in graph theory and computer science 
H . and pose similar hard counting problems. 

\^ ' The chromatic polynomial is also of direct relevance to statistical physics as it is equivalent to the zero-temperature 
partition function of the Potts antiferromagnet [3, 4]: The Potts model [3] constitutes a paradigmatic characterization 
of systems of interacting electromagnetic moments or spins, where each spin can be in one out of g > 2 states; it thus 
Q \ generalizes the Ising model where q = 2. For antiferromagnetic interactions, neighboring spins tend to disalign such 
'— that at zero temperature, the partition function of the Potts antiferromagnet counts the number of ground states 
^ ■ of a spin system just as the chromatic polynomial counts the number of proper colorings of the same graph. For 
^ , sufficiently large q there are many system configurations in which all pairwise interaction energies are minimized at 
' zero temperature. Indeed, these systems exhibit a large number of disordered ground states that is exponentially 
T— I \ increasing with system size. Thus the Potts model exhibits positive ground state entropy, an exception to the 
^) ' third law of thermodynamics. Experimentally, complex disordered ground states and related residual entropy at low 
'nI" \ temperatures have been observed in various systems [5-10]. 

. Although there are several analytical approaches to find chromatic polynomials for families of graphs and to bound 
^) ' their values [1-4, 11-17], there is no closed form solution to this counting problem for general graphs. Algorithmically 
On . it is hard to compute the chromatic polynomial, because the computation time in general increases exponentially 
' with the number of edges in the graph [18]. It also strongly depends on the structure of the graph and rapidly 
increases with the graph's size, and the degrees of its vertices, cf. [19-21]. Therefore, most studies on chromatic 
polynomials up to date have focused on small graphs and famihes of graphs of simple structure and low vertex 
. _ degrees, e.g. two-dimensional lattice graphs [15-17] (an interesting recent attempt to analytically study simple cubic 
■ lattices considered strips with reduced degrees [11]). In fact, it is not at all straightforward to computationally access 
larger graphs with more involved structure, including physically relevant three-dimensional lattice graphs. Finding the 
chromatic polynomial of a graph thus constitutes a challenging, computationally hard problem of statistical physics, 
graph theory and computer science (cf. [18, 19]). 

Below we present a novel, efficient method to compute chromatic polynomials of larger structured graphs. Repre- 
senting a chromatic polynomial as a zero temperature partition function of the Potts antiferromagnet we transform 
the computation into a local and vertex-oriented, ordered pattern matching problem which we then implement using 
appropriate computer algebra. In contrast to conventional top-down methods that represent and process the entire 
graph (and many modified copies thereof) from the very beginning, the new method presented here works through 
the graph bottom-up and thus processes comparatively small local parts of the graph only. 

Consider a graph G that is defined by a set of N vertices i G V = {1, . . . , N} and a set of M :— \E\ edges {i, j} & E, 
each edge joining two vertices i and j which are then called adjacent or neighboring. This graph is said to be (properly) 
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Figure 1: The bottom- up, local and vertex-oriented nature of the algorithm. Sequential processing of edges (dark) adjacent 
to vertices (a) i = 1, (b) i = 2, and (c) i = 3. The remainder graph (light) is not affected when processing vertex i. Partial 
partition functions Zi-i computed so far (before each vertex step i) are shown in square brackets. 



g-colored if every vertex is given one out of q colors {1,2,..., q\ such that every two adjacent vertices have different 
colors. The number of g-colorings of a graph G is expressed by its chromatic polynomial P{G, q), a polynomial in q 
of order N [2]. 

The deletion-contraction theorem of graph theory [2] suggests a simple algorithm to compute the chromatic poly- 
nomial of a given graph recursively. In principle, this algorithm works for arbitrary graphs and is therefore, with 
certain improvements, implemented in general-purpose computer algebra systems such as Mathematica [22-24] and 
Maple [25] (cf. also [26, 27]). However, applying the theorem recursively the chromatic polynomial of exponentially 
many graphs must be found, the (weighted) sum of which yields the chromatic polynomial of the original graph. This 
reflects how hard the problem is algorithmically and severely restricts the applicability of computational methods, in 
particular if they employ standard top-down processing. 

We now describe our novel algorithm. It is based on the antiferromagnetic (J < 0) Potts model [3, 4] with 
Hamiltonian 

giving the total energy of the system in state cr = {ai, . . . , ajy). Here individual spins cri can assume q different values 
CTi G {1, . . . , q}, generalizing he Ising model (g — 2) [28]. Two spins cr^ and aj on the graph G interact if and only if 
they are neighboring, {i,j} £ E, and in the same state, at = aj, i.e. the Kronecker-delta is Sa-io-j = 1 (otherwise, for 
any pair ai ^ aj, it is Saicr, = 0). Thus the total interaction energy is minized if all pairs of neighboring spins are in 
different states. 

The partition function Z{G,q,T) = J2cF'^^Pi~l^^i^)) positive temperature T = (A:b/3)"^, where /cb is the 
Boltzmann constant, can be represented as 

ZiG,q,T)^J2 n (l + ^-^-.-J (2) 
<T {t,j}eE 

where v = exp(/3J) — 1 e (—1,0]. In the limit T ^ (implying (3 J — > —00 and thus v —1) this partition function 
counts the number of ways of arranging the spins a such that no two adjacent spins are in the same state. Thus the 
zero temperature partition function (2) exactly equals [4] the chromatic polynomial 



P{G,q) = lim Z{G,q,T) 



(3) 
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on the same graph G leading to the representation [4, 12] 

P(G,g)= ^ n (l-<5.,.,) (4) 

0-1=1 CTN=l{j,j}eS 

of the chromatic polynomial in terms of sums over products of Kronecker-deltas. 

The algorithm exploits this representation by expanding the products in (4) and symbolically evaluating the right 
hand side vertex by vertex (cf. Fig. 1), considering each individual sum as an operator. This operator interpre- 
tation relies on a recently studied algebraic structure of expressions containing Kronecker-deltas [29]. Here such an 
operator has the simple actions 



= 1, 



(6) 



and for an arbitrary number r S No of factors, 

if the jp, p € {1, . . . , r}, are pairwise distinct and all jp ^ k. 
For illustration consider the chromatic polynomial 

9 9 « 

P(G,g)= E E E(l-'^-i-^)(l-'^-i-3)(l-'5o203) (9) 

(73 = 1 (72 = 1 (71 = 1 

of a triangular (complete) graph comprised of N = 3 vertices and M = 3 edges. We start at vertex i = 1 by expanding 
the relevant product 

Pi = {1 - 6„,„,){1 - S^,^,) (10) 

— 1 ^(71(72 ^(71(73 ~t~ ^(71 (72 *^<7l (73 (-^-^) 

that is comprised of all factors that contain ai. (We note that already BirkhofF [1] in 1912 used closely related 
expansions to theoretically derive an alternative representation of chromatic polynomials.) Symbolically applying the 
above replacement rules (8) yields a "partial partition function" 

Q 

zi = E^i (12) 

(71=1 

= g- 2 + 5^,^3, (13) 

and thus P{G, q) = X]ct3=i Sct2=i ^'^^^ ~ '^<^2<73)- Proceeding with the vertices i = 2 and i = 3 in a similar fashion, we 
obtain p2 = {q — 2 + Sa2a3){^ — ^(72(73); -^2 = 5^^2=1 = {q — 1)(Q' ~ 2), = Z2, and reduce the chromatic polynomial 
to the final result P{G, q) = Z3 ^ 12a3=i P'^ = lil ^ !)('? ^ 2), successively. 

For a general graph G on N vertices, the algorithm is analogous to the example. First define zq = 1. Then, passing 
through the vertices from i = 1 sequentially up to i = N, 

1. construct and expand pi = Zi-i Oii j}e£;(-'- ~ ^<Ti<jj ) where the product is over all edges incident to i that have 
not been considered before, i.e. j > i; 

2. symbolically evaluate the sum Zi = pi applying the simple rules (5)-(8) given above. 
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Figure 2: Computation time t (in ms) for square lattice samples of sizes 2 x n (main panel) and nxn (inset) increases with the 
number of edges M of the graph. The new method (•) drastically outperforms standard methods used in Mathcmatica (O) 
and Maple (□), both with respect to the scaling of the algorithm (quantified by the local slope) and the absolute effectiveness 
(quantified by the absolute times needed), even for moderately sized graphs. 

These operations are local and vertex oriented in the sense that they jointly consider ah edges {i,j} incident to an 
individual vertex i at any one time. A major advantage of this bottom-up algorithm is that all edges that are not 
currently processed are kept outside the computations until they are needed, quite in contrast to standard top-down 
deletion-contraction algorithms. If a graph G has a layered structure, 



with layers H,^, constituting samples of periodic lattices or aperiodic graphs, the vertices i are selected (i.e. numbered) 
layer by layer such that the operations only affect a particularly small portion of the graph at once (Fig. 1). These 
graphs have bounded tree widths, cf. [30, 31]. More generally, vertices are numbered appropriately beforehand, for 
instance, using minimal band width of the graph as a heuristic criterion [32]. 

The computation of the chromatic polynomial has been reduced to a process of alternating expansion of expressions 
and symbolically replacing terms in an appropriate order. In the language of computer science, these operations arc 
represented as the expanding, matching, and sorting of patterns, making the algorithm suitable for computer algebra 
programs optimized for pattern matching. 

To fully exploit the capabilities of this algorithm, we implemented it using the language Form [33, 34] which is 
specialized to large scale symbolic manipulation problems and as such a successful standard tool for, e.g., Fcynman 
diagram evaluation in precision high-energy physics [35, 36]. 

A practicahy relevant measure for the speed of our method is the total CPU time t it needs for a specific calculation. 
For hard counting problem, one generally expects an exponential increase t « Aexp(am) with the size m ^ 1 of the 
problem, here defined as the number rn = M of edges for chromatic polynomials of general graphs. For graphs with 
bounded tree width [30, 31] the solution time of the counting problem typically only grows exponentiahy with the 
width of the graph, i.e. in the square of the number of vertices in the subgraphs H,y. The factor a in the exponent 
determines the scaling of the computational time with problem size and measures the efficiency of the algorithm, 
whereas the prefactor A fixes the absolute time needed and depends, among others, on the software environment and 
hardware used. 

To compare our method to existing ones, we first computed chromatic polynomials of samples of the two-dimensional 
square lattice with free boundary conditions (2 x n strips that have M = 3n — 2 edges and nxn patches that have 
M = 2n(n — 1) edges). The total computation times t have been measured as a function of M for the new method 
as weh as for the standard methods used in Mathematica [22, 23] and Maple [25], respectively. Figure 3 shows that 
the scaling a of the algorithm (given by the local slope of the data points in the logarithmic plot) of our new method 
is markedly better than the one found for the standard deletion-contraction methods. This implies that the new 
method outperforms these standard computational methods in the absolute computation time by several orders of 
magnitude already for moderately sized graphs (e.g. about six orders of magnitude for n = 10 i.e. M = 28). With 
increasing graph size, the advantages of our method become more pronounced. For example, for the 2 x 100 strip of 
the square lattice (M = 298 edges) the pattern matching method needs a computation time of the order of t ~ 2s 
whereas extrapolation of the data shown in Fig. 3 indicates that the same problem is not computationahy accessible 
using the standard deletion-contraction methods implemented in Mathematica. Second, in contrast to transfer matrix 
or other analytical recurrence methods [12, 13, 15], the above method also works in a simple way for graphs with 
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Figure 3: Computation time t (in ms) for randomly edge-diluted square lattice samples with next nearest neighbor (nnn) 
interactions (displayed as a cartoon inset in panel (a)). The new method (filled symbols) strongly outperforms standard 
deletion-contraction method (Mathematica, open symbols), (a) Computation time vs. the actual total number M of edges for 
samples of 3 x n vertices where each original edge has been deleted independently with probability p = 0.1 (circles) p = 0.2 
(triangles) and p = 0.5 (squares), (b) Computation time vs. the fraction / of edges present. Edges are sequentially randomly 
removed independently, starting from an undiluted 3x5 square lattice graph (/ = 1) with nnn interactions; the graph is 
diluted until before it becomes disconnected below / « 0.4. Whereas the computation times using the new method are still at 
fluctuation level (t < 10~^s), standard methods take at factor of 10^ to lO'^ longer. 



non-identical subgraphs H,^, such as randomly diluted lattices. The same comparison for randomly diluted 3 x n 
square lattice samples with next-nearest neighbor interactions (Fig. 3) confirms the pronounced outperformance and 
moreover illustrates the general applicability of our method, also compared to recursive analytical methods. 

As a first application to an open hard counting problem, we now turn to three-dimensional lattices of direct physical 
relevance. First, we consider nxnx n samples of the simple cubic lattice with free boundary conditions, which have 
N = vertices and M = 3n^{n — 1) edges. We found chromatic polynomials up to n = 4 (A'' = 64, M = 144). A 
representation of the chromatic polynomial P{G, q) in terms of its complex zeroes qi, . . . ,qN is shown in Fig. 4a 
for n = 3 and n = 4. We further consider simple cubic lattice strips that extend in the diagonal (111) direction with 
periodic boundaries in the two other (transverse) directions. This keeps the number of vertices within one layer low at 
the same time allowing for a large number Nc of vertices with the same degree (equal to six) as vertices in the infinite 
lattice, a fact that is heuristically known to be essential for a rapid convergence towards the thermodynamic limit 
{N,M — > oo). The largest three-dimensional sample graphs shown in Figure 4b have N = 384 vertices, M = 1128 
edges and a fraction Nc/N = 368/384 sa 0.96 of vertices with correct degree, (as compared to Nc/N = 8/64 w 0.13 
for the 4x4x4 sample extending along the Cartesian axes and to A'c = for previous attempts to address three- 
dimensional lattices). The computation time was approximately 11 hours on a single Linux machine with an Intel 
Pentium 4, 2.8 GHz-32 bit processor. 

In summary we have presented a novel method to calculate chromatic polynomials of graphs. Using the partition 
function representation, it proceeds vertex by vertex employing an a priori reduction to local operations only and 
is thus particularly suited for graphs exhibiting a layered structure. The method combines a symbolic bottom-up 
algorithm, which is based on systematic term-wise expansion and pattern matching, with an appropriate computer 
algebra program [33, 34]. Our method is applicable to general types of graphs, including graphs with bounded 
and unbounded tree widths as well as randomized graphs. We demonstrated by several sets of examples that it 
drastically outperforms existing standard methods for all these types of graphs. As a practical apphcation, we 
computed chromatic polynomials for samples of the simple cubic lattice, for the first time computationally accessing 
three-dimensional lattices of physical relevance. 

Since the main ideas underlying our method are simple to apply, they may be generalized in a straightforward 
way and also be transferred to other chahenging counting problems. Among others, one may compute quantitative 
measures relevant in computer science that give information about the communication capabilities of a network, 
such as (i) the flow polynomial and (ii) the reliability polynomial [41, 42]. It is equally possible to determine (iii) 
ferromagnetic and (iv) positive temperature partition functions of statistical physics [38-40] and, equivalently, (v) the 
Tutte polynomial (of two variables q and v, cf. Eq. 2), valuations of which directly result in the number of spanning 
subgraphs, the number of spanning trees, and other invariants of a graph [12, 13]. Of course, applications in graph 
theory may include studies of families of graphs where computational results seemed impossible so far, because the 
computational effort is substantially reduced. As the new method yields exact results not only for the final solution 
(in our examples, the chromatic polynomial) but also in the intermediate steps (the partial partition functions above), 
it may moreover be combined with analytical tools [11, 15-17, 38, 40, 42] to obtain unprecedented results for various 
classes of graphs. Finally, the method can easily be implemented in parallel computations. Taken together, the 
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Figure 4: Complex zeroes representing chromatic polynomials of samples of the simple cubic lattice (a) n x n x n Cartesian 
samples with free boundary conditions, and (b) 2 x 4 x na diagonal samples with periodic transverse boundary conditions, with 
up to M = 1128 edges. 

novel bottom-up pattern matching algorithm combined with specialized computer algebra presented here constitutes 
a promising starting point to access a number of challenging, computationally hard counting problems from statistical 
physics, graph theory and computer science. 
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